﻿import numpy as np
import math
import matplotlib.pyplot as plt
from mpl_toolkits import mplot3d
from stl import mesh
import pyopencl as cl
import pylas
import time

"""
Описывает основной функционал работы модели канала ЛЛС
"""


class Lidar:
    def __init__(self, fov_h, fov_v, step_h, step_v):
        start_time = time.time()
        """
        :param fov_h: ширина сканирования по азимуту
        :param fov_v: ширина сканирования по углу места
        :param step_h: разрешение по азимуиу
        :param step_v: разрешение по углу места
        """
        self.fov_h = fov_h
        self.fov_v = fov_v

        self.step_h = step_h
        self.step_v = step_v

        self.quantity_rays_v = np.floor(self.fov_v / self.step_v).astype(np.int32)
        self.quantity_rays_h = np.floor(self.fov_h / self.step_h).astype(np.int32)

        self.dir_data = np.zeros((self.quantity_rays_v, self.quantity_rays_h, 3), dtype=np.float32)
        self.fov_v = self.quantity_rays_v * self.step_v
        self.fov_h = self.quantity_rays_h * self.step_h

        self.central_dir = np.array([0, 1, 0], dtype=np.float32)  # вектор нормали окна
        self.position = np.array([0, 0, 0], dtype=np.float32)
        self.lu_dir = np.array(
            [self.central_dir[1] * math.tan(math.radians(-self.fov_h / 2)), self.central_dir[1],
             self.central_dir[1] * math.tan(math.radians(self.fov_v / 2))], dtype=np.float32)  # левый верхний угол матрицы
        self.ld_dir = np.array(
            [self.central_dir[1] * math.tan(math.radians(-self.fov_h / 2)), self.central_dir[1],
             self.central_dir[1] * math.tan(math.radians(-self.fov_v / 2))], dtype=np.float32)  # левый нижний угол матрицы
        self.ru_dir = np.array(
            [self.central_dir[1] * math.tan(math.radians(self.fov_h / 2)), self.central_dir[1],
             self.central_dir[1] * math.tan(math.radians(self.fov_v / 2))], dtype=np.float32)  # правый верхний угол матрицы
        self.rd_dir = np.array(
            [self.central_dir[1] * math.tan(math.radians(self.fov_h / 2)), self.central_dir[1],
             self.central_dir[1] * math.tan(math.radians(-self.fov_v / 2))], dtype=np.float32)  # правый нижний угол матрицы
        # print((self.ru_dir-self.lu_dir)/2)

        self.dir_data[:] = self.ld_dir
        self.width_win = self.ru_dir[0] - self.lu_dir[0]  # ширина матрицы
        self.height_win = self.ru_dir[2] - self.rd_dir[2]  # высота матрицы
        self.single_qr_v = 20
        self.single_qr_h = 20
        self.step_win_v = self.height_win / self.quantity_rays_v
        self.step_win_h = self.width_win / self.quantity_rays_h

        # формируем матрицу лучей
        h_step = (np.linspace(self.step_win_h / 2, self.width_win + self.step_win_h / 2, self.quantity_rays_h))
        v_step = (np.linspace(self.step_win_v / 2, self.height_win + self.step_win_v / 2, self.quantity_rays_v))
        h_step = h_step[np.newaxis, :]
        v_step = v_step[:, np.newaxis]
        self.dir_data[:, :, 0] += h_step
        self.dir_data[:, :, 2] += v_step
        self.dir_data = self.dir_data.reshape((self.quantity_rays_v * self.quantity_rays_h, 3))
        self.dir_matrix = np.flipud(np.reshape(self.dir_data, (self.quantity_rays_v, self.quantity_rays_h, 3)))
        self.cloud_points = np.array([])
        self.dp = np.array([])
        self.single_dp = np.array([])
        self.matrix_delay = np.array([])
        self.ih = np.array([])
        self.ih_time = np.array([])
        self.out_signal = np.array([])
        self.out_signal_time = np.array([])
        self.t_d = .0
        print('Time init', time.time() - start_time)

    def rotate_window(self, axis, angle):
        if axis == 'x':
            rotate_matrix = np.array(
                [[1, 0, 0],
                 [0, math.cos(math.radians(angle)), -math.sin(math.radians(angle))],
                 [0, math.sin(math.radians(angle)), math.cos(math.radians(angle))]])
        elif axis == 'z':
            rotate_matrix = np.array(
                [[np.cos(np.deg2rad(angle)), -np.sin(np.deg2rad(angle)), 0],
                 [np.sin(np.deg2rad(angle)), np.cos(np.deg2rad(angle)), 0],
                 [0, 0, 1]])
        elif axis == 'y':
            rotate_matrix = np.array(
                [[math.cos(math.radians(angle)), 0, math.sin(math.radians(angle))],
                 [0, 1, 0],
                 [-math.sin(math.radians(angle)), 0, math.cos(math.radians(angle))]])
        else:
            print('Проверьте ось вращения')
            rotate_matrix = np.array([[0, 0, 0],
                                      [0, 0, 0],
                                      [0, 0, 0]])
        self.dir_data[:] = np.dot(self.dir_data[:], rotate_matrix.transpose())

    def set_position(self, position):
        """
        Переместить ЛЛС в точку с координатами x,y,z
        :param x:
        :param y:
        :param z:
        :return:
        """
        position = position.astype(np.float32)
        self.position = position

    def move_position(self, delta):
        """
        Переместить ЛЛС на delta
        :param delta:
        :return:
        """
        self.position += delta

    def show_matrix(self, path):
        figure = plt.figure()
        scene = mesh.Mesh.from_file(path)
        axes = mplot3d.Axes3D(figure)
        axes.scatter(self.central_dir[0], self.central_dir[1], self.central_dir[2], c='green')
        # axes.scatter(self.lu_dir[0], self.lu_dir[1], self.lu_dir[2], c='orange')
        # axes.scatter(self.ru_dir[0], self.ru_dir[1], self.ru_dir[2], c='orange')
        # axes.scatter(self.ld_dir[0], self.ld_dir[1], self.ld_dir[2], c='orange')
        # axes.scatter(self.rd_dir[0], self.rd_dir[1], self.rd_dir[2], c='orange')

        axes.scatter(self.dir_data[1][0], self.dir_data[1][1], self.dir_data[1][2])
        axes.scatter(self.dir_data[-1][0], self.dir_data[-1][1], self.dir_data[-1][2])
        axes.scatter(self.position[0], self.position[1], self.position[2], c='black')
        scale = scene.points.flatten()
        axes.auto_scale_xyz(scale, scale, scale)
        plt.show()

    def create_plot_dp(self, image_array, mode=None, figure=None, canvas=None):
        """
        Построение дальностной матрицы
        :param image_array: дальностная матрица
        :param mode: 'single' для одного канала
        :param figure:
        :param canvas:
        :return:
        """
        if mode == 'single':
            image_array = np.reshape(image_array, (self.single_qr_v, self.single_qr_h))
        else:
            image_array = np.reshape(image_array, (self.quantity_rays_v, self.quantity_rays_h))
        image_array = np.flipud(image_array)


        # np.savez('./result_dp/1.npz', image_array) # save data manually


        if figure and canvas:
            figure.clear()
            ax = figure.add_subplot(111)
            plot = ax.imshow(image_array)
            figure.colorbar(plot)
            canvas.draw()
        else:
            plt.figure()
            plt.imshow(image_array)
            plt.colorbar()
            plt.show()

    def scan(self, polygons, mode=None, directions=None):
        """
        Определение пересечения лучей и сцены, формирование облака точек и дальностной матрицы
        :param polygons: полигоны сцены
        :param mode: 'single' для одного канала
        :param directions: массив лучей, указывается в случае сканирования одного канала
        :return:
        """
        self.cloud_points = np.array([])
        polygons = polygons.astype(np.float32)
        n_poly = len(polygons)
        ctx = cl.create_some_context()
        queue = cl.CommandQueue(ctx)
        mf = cl.mem_flags
        src = """
        float ray_triangle_intersection(float3 v0, float3 v1, float3 v2, float3 orig, float3 dir)
        {
            float3 e1 = v1 - v0;
            float3 e2 = v2 - v0;
            float3 pvec = cross(dir, e2);
            float det = dot(e1, pvec);
            //float eps = 0.0000001;
            if (det < 1e-8f && det > -1e-8f) {
                return 0.;
            }
            float inv_det = 1 / det;
            float3 tvec = orig - v0;
            float u = dot(tvec, pvec) * inv_det;
            if (u < 0 || u > 1) {
                return 0.;
            }
            float3 qvec = cross(tvec, e1);
            float v = dot(dir, qvec) * inv_det;
            if (v < 0 || u + v > 1) {
                return 0.;
            }
            float t = dot(e2, qvec) * inv_det;
            return t;
        }
        
         __kernel void rti(
         __global const float *start_cl, __global const float *directions_cl, __global const float *poly_cl,
         int const n_poly, __global float *res)
         {
          int gid = get_global_id(0);
          float3 orig = (float3)(start_cl[0], start_cl[1], start_cl[2]);
          float3 dir = (float3)(directions_cl[gid*3], directions_cl[gid*3+1], directions_cl[gid*3+2]);
          float t_min = 1000;
          
          for (int poly_i = 0; poly_i < n_poly; poly_i++){
            float3 v0 = (float3)(poly_cl[poly_i*9], poly_cl[poly_i*9+1], poly_cl[poly_i*9+2]);
            float3 v1 = (float3)(poly_cl[poly_i*9+3], poly_cl[poly_i*9+4], poly_cl[poly_i*9+5]);
            float3 v2 = (float3)(poly_cl[poly_i*9+6], poly_cl[poly_i*9+7], poly_cl[poly_i*9+8]);
            float t = ray_triangle_intersection(v0, v1, v2, orig, dir);
            if (t < t_min && t!=0 && t>1e-8f){
                t_min = t;
            }  
          }
          
          if (t_min != 1000){
            float3 point = orig+t_min*dir;
            res[gid*4] = point.x;
            res[gid*4+1] = point.y;
            res[gid*4+2] = point.z;
            res[gid*4+3] = length(point-orig);
          }
          else if (t_min == 1000){
            res[gid*4] = 0;
            res[gid*4+1] = 0;
            res[gid*4+2] = 0;
            res[gid*4+3] = 0;
          }  
         }
        """
        prg = cl.Program(ctx, src).build()
        start = self.position.copy()
        # print(start.dtype)
        start_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=start)
        if mode == 'single':
            direction_data = directions
        else:
            direction_data = self.dir_data
        direction_data_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=direction_data)
        poly_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=polygons)
        res_g = cl.Buffer(ctx, mf.WRITE_ONLY, direction_data.nbytes * 2)
        knl = prg.rti  # Use this Kernel object for repeated calls
        knl(queue, (len(direction_data), 1), None, start_cl, direction_data_cl, poly_cl, np.int32(n_poly), res_g)
        size_res = (np.empty([len(direction_data), 4])).astype(np.float32)
        # if mode == 'single':
        # self.single_dp = np.empty_like(size_res_dis)
        # cl.enqueue_copy(queue, self.single_dp, dis_array)
        # self.single_dp[self.single_dp == 0] = np.inf
        # else:
        #     self.dp = np.empty_like(size_res_dis)
        #     cl.enqueue_copy(queue, self.dp, dis_array)
        #     self.dp[self.dp == 0] = np.inf
        res_np = np.empty_like(size_res)
        cl.enqueue_copy(queue, res_np, res_g)
        # cl.enqueue_copy(queue, self.dp, dis_array)
        self.dp = res_np[:, 3].copy()
        self.dp[self.dp == 0] = np.inf
        cloud_points = []
        for point in res_np:
            if (point != 0.).any():
                cloud_points.append(point)
        self.cloud_points = np.asarray(cloud_points)

    def scan_bvh(self, aabblist, triangleslist, tristartlist, ntrilist, childlist, mode=None, directions=None):
        self.cloud_points = np.array([])
        start = self.position.copy()
        ctx = cl.create_some_context()
        queue = cl.CommandQueue(ctx)
        mf = cl.mem_flags
        src = """          
                
                int isempty(int *top, int *MAXSIZE, int *stack) {                
                   if(*top == -1)
                      return 1;
                   else
                      return 0.0f;
                }
                   
                int isfull(int *top, int *MAXSIZE, int *stack) {                
                   if(*top == *MAXSIZE)
                      return 1;
                   else
                      return 0.0f;
                }
                
                int peek(int *top, int *MAXSIZE, int *stack) {
                   return stack[*top];
                }
                
                int pop(int *top, int *MAXSIZE, int *stack) {
                   int data;
                   if(!isempty(top, MAXSIZE, stack)) {
                      data = stack[*top];
                      *top = *top - 1;   
                      return data;
                   }
                   else{
                   return data;
                   }
                }
                
                void push(int data, int *top, int *MAXSIZE, int *stack) {
                   if(!isfull(top, MAXSIZE, stack)) {
                      *top = *top + 1;   
                      stack[*top] = data;
                   }
                }
                
                
                float ray_triangle_intersection(float3 v0, float3 v1, float3 v2, float3 orig, float3 dir)
                {
                    float3 e1 = v1 - v0;
                    float3 e2 = v2 - v0;
                    float3 pvec = cross(dir, e2);
                    float det = dot(e1, pvec);
                    //float eps = 0.0000001;
                    if (det < 1e-8f && det > -1e-8f) {
                        return 0.0f;
                    }
                    float inv_det = 1 / det;
                    float3 tvec = orig - v0;
                    float u = dot(tvec, pvec) * inv_det;
                    if (u < 0 || u > 1) {
                        return 0.0f;
                    }
                    float3 qvec = cross(tvec, e1);
                    float v = dot(dir, qvec) * inv_det;
                    if (v < 0 || u + v > 1) {
                        return 0.0f;
                    }
                    float t = dot(e2, qvec) * inv_det;
                    return t;
                }

                bool ray_box_intersection(float3 orig, float3 dir, float3 boxMin, float3 boxMax)
                {
                    float3 inv_dir = 1/dir;
                    float lo = inv_dir.x*(boxMin.x - orig.x);
                    float hi = inv_dir.x*(boxMax.x - orig.x);
                    float tmin  = min(lo, hi);
                    float tmax = max(lo, hi);
                    float lo1 = inv_dir.y*(boxMin.y - orig.y);
                    float hi1 = inv_dir.y*(boxMax.y - orig.y);
                    tmin  = max(tmin, min(lo1, hi1));
                    tmax = min(tmax, max(lo1, hi1));
                    float lo2 = inv_dir.z*(boxMin.z - orig.z);
                    float hi2 = inv_dir.z*(boxMax.z - orig.z);
                    tmin  = max(tmin, min(lo2, hi2));
                    tmax = min(tmax, max(lo2, hi2));
                    return (tmin <= tmax) && (tmax > 0.f);
                }
                

                 __kernel void rti(
                 __global const float *start_cl, __global const float *directions_cl, __global const float *poly_cl,
                 __global const float *aabblist_cl, __global const int *tristartlist_cl,
                 __global const int *ntrilist_cl, __global const int *childlist_cl,
                 __global float *res)
                 {
                  int gid = get_global_id(0);
                  int MAXSIZE = 100;       
                  int stack[100];     
                  int top = -1;
                  float3 orig = (float3)(start_cl[0], start_cl[1], start_cl[2]);
                  float3 dir = (float3)(directions_cl[gid*3], directions_cl[gid*3+1], directions_cl[gid*3+2]);
                  float tmin = 10000; 
                  push(0, &top, &MAXSIZE, stack);
                  while (isempty(&top, &MAXSIZE, stack) != 1){
                    int n_node = pop(&top, &MAXSIZE, stack);
                    float3 boxMin = (float3)(aabblist_cl[n_node*6], aabblist_cl[n_node*6+1], aabblist_cl[n_node*6+2]);
                    float3 boxMax = (float3)(aabblist_cl[n_node*6+3], aabblist_cl[n_node*6+4], aabblist_cl[n_node*6+5]);
                    int child1 = childlist_cl[n_node*2];
                    int child2 = childlist_cl[n_node*2+1];
                    bool rbi = ray_box_intersection(orig, dir, boxMin, boxMax);
                    if (rbi == true && child1 != -1 && child2 != -1){
                        push(child2, &top, &MAXSIZE, stack);
                        push(child1, &top, &MAXSIZE, stack);
                    }
                    else if (rbi == true && child1 == -1 && child2 == -1){
                        int tristart = tristartlist_cl[n_node];
                        int ntri = ntrilist_cl[n_node];
                        for (int i = tristart; i < tristart+ntri; i++){
                            float3 v0 = (float3)(poly_cl[i*9], poly_cl[i*9+1], poly_cl[i*9+2]);
                            float3 v1 = (float3)(poly_cl[i*9+3], poly_cl[i*9+4], poly_cl[i*9+5]);
                            float3 v2 = (float3)(poly_cl[i*9+6], poly_cl[i*9+7], poly_cl[i*9+8]);
                            float t = ray_triangle_intersection(v0, v1, v2, orig, dir);
                            if (t < tmin && t!=0 && t>1e-8f){
                                tmin = t;

                            }  
                        }
                    }
                    if (tmin!=10000){
                        float3 point = orig+tmin*dir;
                        res[gid*4+3] = length(point-orig);
                        //res[gid*4+3] = tmin;
                        res[gid*4] = point.x;
                        res[gid*4+1] = point.y;
                        res[gid*4+2] = point.z;
                    }
                    else if(tmin == 10000){
                        
                        res[gid*4+3] = 0;
                        res[gid*4] = 0;
                        res[gid*4+1] = 0;
                        res[gid*4+2] = 0;                    
                    }   
                  }
                 }
                """
        prg = cl.Program(ctx, src).build()
        aabblist_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=aabblist)
        tristartlist_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=tristartlist)
        ntrilist_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=ntrilist)
        childlist_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=childlist)
        start_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=start)
        if mode == 'single':
            direction_data_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.USE_HOST_PTR, hostbuf=directions)
            res_np = (np.empty([len(directions), 4])).astype(np.float32)
            res_g = cl.Buffer(ctx, mf.WRITE_ONLY, directions.nbytes * 2)
        else:
            direction_data_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.USE_HOST_PTR, hostbuf=self.dir_data)
            res_np = (np.empty([len(self.dir_data), 4])).astype(np.float32)
            res_g = cl.Buffer(ctx, mf.WRITE_ONLY, self.dir_data.nbytes * 2)

        poly_cl = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=triangleslist)
        knl = prg.rti
        knl(queue, (len(self.dir_data), 1), None, start_cl, direction_data_cl, poly_cl, aabblist_cl, tristartlist_cl,
            ntrilist_cl, childlist_cl, res_g)
        cl.enqueue_copy(queue, res_np, res_g)
        if mode == 'single':
            self.single_dp = res_np[:, 3].copy()
            self.single_dp[self.single_dp == 0] = np.inf
        else:
            self.dp = res_np[:, 3].copy()
            self.dp[self.dp == 0] = np.inf
        cloud_points = []
        for point in res_np:
            if (point != 0.).any():
                cloud_points.append(point)
        self.cloud_points = np.asarray(cloud_points)

    def write_las(self, name):
        """
        Запись облака точек в файл
        :param name: имя файла
        :return:
        """
        las = pylas.create()
        las.x = self.cloud_points[:, 0]
        las.y = self.cloud_points[:, 1]
        las.z = self.cloud_points[:, 2]
        las.write(name)

    def get_ih(self, dp):
        """
        Получение импульсной характеристики
        :param dp: дальностная матрица
        :return:
        """
        matrix_delay = 2 * dp / 3e8
        delay_flatten = matrix_delay.flatten()
        delay_flatten = delay_flatten[delay_flatten != np.inf]
        delay_flatten = np.append(delay_flatten, [0])
        delay_flatten = np.append(delay_flatten, [5e-6])
        delay_flatten.sort()
        # print('Матрица задержек', delay_flatten)
        # print('Time discrete', self.t_d)

        self.ih, bins = np.histogram(delay_flatten, int((np.amax(delay_flatten) - np.amin(delay_flatten)) / self.t_d))
        self.ih_time = np.empty_like(self.ih).astype(np.float32)
        for i in range(len(self.ih)):
            self.ih_time[i] = (bins[i] + bins[i + 1]) / 2
        self.ih = self.ih / (self.single_qr_v * self.single_qr_h)

    def create_plot_ih(self, figure=None, canvas=None):
        """
        Построение графика ИХ
        :param figure:
        :param canvas:
        :return:
        """
        if figure and canvas:
            figure.clear()
            ax = figure.add_subplot(111)
            ax.plot(self.ih_time, self.ih)
            ax.grid()
            ax.set_xlabel('Время, с')
            ax.set_ylabel('Амплитуда')
            ax.set_title('Импульсная характеристика')
            canvas.draw()

    def create_plot_out_signal(self, in_signal, figure=None, canvas=None):
        """
        Формирование выходного сигнала и построение его графика
        :param in_signal: входной сигнал
        :param figure:
        :param canvas:
        :return:
        """
        if figure and canvas:
            figure.clear()
            ax = figure.add_subplot(111)
            noise = np.random.normal(0, 10, len(self.ih) + len(in_signal) - 1)
            self.out_signal = np.convolve(self.ih, in_signal) + noise
            # self.out_signal = np.convolve(self.ih, in_signal+np.random.poisson(1, len(in_signal)))
            self.out_signal_time = np.linspace(0, (len(self.ih) + len(in_signal) - 1) * self.t_d,
                                               len(self.ih) + len(in_signal) - 1)
            ax.plot(self.out_signal_time, self.out_signal)
            # ax.plot(self.out_signal + np.random.normal(0, 10, np.size(self.out_signal)))
            ax.grid()
            ax.set_xlabel('Время, с')
            ax.set_ylabel('Амплитуда')
            ax.set_title('Выходной сигнал')
            canvas.draw()
        else:
            noise = np.random.normal(0, 0.01, len(self.ih) + len(in_signal) - 1)
            # noise = np.random.poisson(0.1, len(self.ih)+len(in_signal)-1)
            self.out_signal = np.convolve(self.ih, in_signal) + noise
            self.out_signal_time = np.linspace(0, (len(self.ih) + len(in_signal) - 1) * self.t_d,
                                               len(self.ih) + len(in_signal) - 1)

    def get_single_dir(self, ch_h, ch_v):
        """
        Формирование массива лучей для одного канала
        :param ch_h: номер канала по азимуту
        :param ch_v: номер канала по углу места
        :return:
        """
        single_central_dir = self.dir_matrix[ch_v][ch_h]
        single_ld_dir = np.array(
            [single_central_dir[0] - self.step_win_h / 2, single_central_dir[1], single_central_dir[2]
             - self.step_win_v / 2])
        single_dir_data = []
        for v in np.linspace(0, self.step_win_v, self.single_qr_v):
            for h in np.linspace(0, self.step_win_h, self.single_qr_h):
                single_dir_data.append(np.array([single_ld_dir[0] + h, single_ld_dir[1], single_ld_dir[2] + v]))
        single_dir_data = np.asarray(single_dir_data).astype(np.float32)
        return single_dir_data

    def get_all_out_signal(self, aabblist, triangleslist, tristartlist, ntrilist, childlist, in_signal):
        """
        формирование выходного сигнала для всех каналов ЛЛС
        :param childlist:
        :param ntrilist:
        :param tristartlist:
        :param triangleslist:
        :param aabblist:
        :param in_signal:
        :return:
        """
        # len_out = len(self.ih) + len(in_signal) - 1
        all_out = np.zeros((199, self.quantity_rays_v, self.quantity_rays_h))
        # print(all_out.shape)
        for v in range(self.quantity_rays_v):
            for h in range(self.quantity_rays_h):
                directions = self.get_single_dir(h, v)
                # print('Матрица направлений сформирована')
                self.scan_bvh(aabblist, triangleslist, tristartlist, ntrilist, childlist, mode='single',
                              directions=directions)
                # print('Сканирование завершено')
                self.get_ih(self.single_dp)
                # print("Импульсная характеристика сформирована")
                self.create_plot_out_signal(in_signal)
                # print("Выходной сигнал рассчитан")
                # print(self.out_signal.shape)
                all_out[:, v, h] = self.out_signal
                # print(self.out_signal)
                print(f'Канал {v}, {h} обработан')
        return all_out
